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flocking behaviour for an experimental interaction rule 
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Computational models of collective behavior in birds has allowed us to infer interaction rules di- 
rectly from experimental data. Using a generic form of these rules we explore the collective behavior 
and emergent dynamics of a simulated swarm. For a wide range of flock size and interaction extent 
(the fixed number of neighbors with which an individual will interact) we find that the computa- 
tional collective is inherently stable — individuals are attracted to one another and will position 
themselves a preferred distance from their fixed neighbors within a rigid lattice. Nonetheless, the 
irregular overall shape of the flock, coupled with the need for individuals on the boundary to move 
towards their neighbors creates a torque which leads the flock to rotate and then meander. We argue 
that this "rolling meander" is a very good proxy for real collective behavior in animal species and 
yet arises from a simple homogeneous and deterministic rule for interaction. Rather than then in- 
troduce leaders — which has already been shown, quite straightforwardly, to drive collective swarms 
such as this — we introduce a small number of "followers". Each follower is bound to consider a 
random fixed individual to be among their neighbors, irrespective of actual metric distance between 
them. We find that the introduction of a small number of such followers causes a phase transition 
that quickly leads to instability in the flock structure (as no stable configuration arises) and the 
previously rigid crystalline interaction among neighbors now becomes fluid: the distance between 
neighbors decreases, the flock ceases to rotate and meanders less. 



Collective behavior occurs widely in nature, and as a 
consequence efforts to understand it have been vigorous 
PQ. Physicists and computer scientists have had their 
own unique interest in this problem as simple mathe- 
matical models for inter-particle interaction have been 
shown to lead to interesting collective behaviors [2] and 
even phase transitions as a function of particle density 
[3]. While biologists approach this problem with efforts 
to build intricate models of inter-individual interactions 
that capture the known biology, physicists have focussed 
on these ad hoc models with simple mechanics that pro- 
duce "naturalistic" behaviors. Recent advances in track- 
ing and monitoring technology now allows for a third 
approach. By observing real flocks it is possible to build 
models of inter-individual interaction directly from the 
data [4] . In [5] we did this for data from a group of hom- 
ing pigeons in free flight. One of the observations of this 
work was that the hierarchical leadership relationship de- 
tected originally [4] was actually supplemented by a more 
egalitarian reciprocity [6]. 

In this communication we take these models built from 
pigeon flight data and ask what features can we observe 
in the models which are necessary to explain observed 
naturalistic flight? And, what dynamical phenomena do 
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these models exhibit? This is a different problem from 
that treated by Vicsek [3 j and Reynolds [2 — for whom 
the model was created with a combination of intuition 
and foresight. It is also different from what has been 
done in biology [H [71 [8] , as we do not ground our model 
directly on theory but rather on experimental observa- 
tion and data. Our work diverges still further from the 
approach taken in biology as we distill from the compu- 
tational model a simple interaction rule consistent with 
both the observed interaction data and collective behav- 
ior. 

The details of our modeling procedure have been de- 
scribed elsewhere [5 and, while foundational, are not di- 
rectly relevant to the current discussion. Nonetheless, 
we build a model of particle-particle interaction where 
the position x[^ of particle i is updated as a function 
/(•) of its current position and velocity and those of 
its m nearest neighbors. The function / is chosen, via 
a numerical procedure [9HTT]. to fit the observed data. 
The exact details of this procedure are not relevant here 
and any of the many alternatives would do just as well. 
Nonetheless, from this model / we can now estimate the 
strength of pairwise interaction between two particles. If 
the distance between two particles i and j is given by 
dij = \\x^ — x^\\ then we define the force of attraction 
h(dij). The form of hj-) can be estimated from the data 
and is shown in Fig. [11(a). 

The exact (that is, exact representation of the model 
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FIG. 1. Pairwise attraction repulsion. [Colour Online] 
Panel (a) depicts the pairwise force of attraction or repulsion 
between two particles based on simulations of the computa- 
tional model built from the data — following methods de- 
scribed extensively in [5]. Panel (b) is the simplified proxy 
for this curve h(x) = — x h e~ x (for a — 5 and b — 2) 
which we will employ in the current work. 



/) interaction rule h is computationally complicated. In 
Fig. [IJb) we depict the proxy we will use here instead: 
h. In addition to being simpler to compute, h has the 
advantage that it is more easily replicable. While we 
could use the exact curve depicted in Fig. []Ja), this is 
not actually necessary. It is only the general shape which 
we wish to capture and hence we use the function 



h[x) = 



x + a 



(1) 



as a computationally advantageous expediency. With 
h(x) providing an attractive/repuslve force for the inter- 
action between two particles based on distance, we are 
now able to state the swarm dynamic which we study 
here. 

(i) 

Let the swarm consist of TV particles such that x t is 
the position and the velocity of particle i at time t. 
Each particle updates its velocity and position based on 
interaction with its k nearest neighbors. Denote by x 1 ^ 3 
the j — th nearest neighbor of particle i at time t. We 
then compute two competing desired velocity vectors for 
each particle 



dij 
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where, for each particle, vy is the desired velocity from 
velocity alignment amount the neighbors and vh is the 
desired velocity based on the pairwise attraction or re- 
pulsion with each neighbor. The actual updated velocity 
is a convex combination of these two desired quantities, 
and the position is then computed from this: 
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where < A < 1 and e is an appropriate integration time 
step. In what follows we have set A = 0.3, but this choice 
is not critical, extensive repeated simulations with other 
intermediate (i.e. A ^ 0, 1) values provide equivalent 
results. 

The qualitative behaviors of simulated swarms follow- 
ing these rules are striking, and have been briefly dis- 
cussed in the abstract. For all values of N (up to 5000) 
and k (1 < k < 28) our simulations indicate that the 
swarm will find a stable configuration — initially ran- 
domly distributed particles (both in space and velocity) 
quickly rearrange themselves into a stable configuration 
reminiscent (for smaller k) of a sub-optimal sphere pack- 
ing (each particle strives to reach its preferred distance 
from the requisite neighbors). The position of each par- 
ticle relative to its neighbors then remains fixed within 
a single cloud. However, since particle on the boundary 
of that cloud have neighbor to only one side of them, 
these particles attempt to move towards their neighbors 
in such a way that the irregular shape of the boundary 
creates a non-zero rotational torque. The cloud rotates 
and then meanders. 

Fascinating as this rigid form of self organization is, we 
prefer to pursue a more quantitative analysis here. To do 
so we define the following six quantities: 

• local separation: the average over all particles in 
the flock of the distance between that point and its 
k nearest neighbors; 

• global separation: the maximum distance between 
any two particles in the swarm; 

• neighbor stability: the average over all particles 
of number of different particles which become the 
nearest neighbor of a given particle over a fixed 
time interval; 

• vorticity: the mean rotation, between successive 
time intervals of each particle around the centre 
of mass of the swarm; 

• travel: the path integral of the centre of mass over 
a time interval; and, 

• movement: the straight line distance travelled over 
the same time period. 

For a wide range of values of k and N we simulated the 
behaviour of swarms of particles from random initial con- 
ditions — random position and velocity. In each case we 
simulated for 20000 time steps and computed the pa- 
rameters defined above over the last 500 (sampling only 
everything 10-th time step) configurations. In Fig. [2] we 
depict the behaviour of each of these six parameters. 

Unsurprisingly, local separation depends primarily on 
the number k (imagine a single ball trying to be equally 
close to k other balls as k increases) and global separa- 
tion depends on N. Nonetheless, we do see for a fixed 
k a slight variance in local separation as function of TV 
— separation is greatest around N = 400. Similarly, 
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(b) global sep. 
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(d) vorticity 
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FIG. 2. Effect of swarm size and number of neighbors 

[Colour Online] Panel (a-f) depicts heat maps and/or contour 
plots for each of the size parameters defined in the text (lo- 
cal separation, global separation, neighbor stability, vorticity, 
travel and movement) for k = 1, 2, ... 28 and 100 < N < 5000. 
Darker coloring corresponds to smaller numerical values. 




FIG. 3. Phase transition [Colour Online] (a) local sep- 
aration, (b) global separation, (c) neighbor stability, and 
(d) vorticity as a function of k = 1,2,... 28 for N = 
80, 90, 100, 200, 300, 400, 500. Behaviour of local and global 
separation is as describe above for Fig. [2] For neighbor stabil- 
ity there is a clear change in stability as k increases the swarm 
ceases to have a rigid equilibrium and the nearest neighbors 
constantly switch, for larger still k(k > 20) the stable rigid 
structure returns. The width of the constant window increases 
with N (over the range of N plotted here): N = 300, 400, 500 
show the largest deviation from a stable rigid structure. Vor- 
ticity shows a clear increase are k exceeds about 6 — beyond 
this point the swarm will start to rotate and continue to do 
so. Significant variability between simulations is particularly 
evident here. 



for a fixed N global separation is decreased by increas- 
ing k: as a larger k introduces a stronger cohesive force: 
nonetheless, there is a slight peak in global separation for 
intermediate values of k (around k = 10, see Fig. [3](b)). 
The plot of neighbor stability is complementary to these 
findings as we see that for small and more evidently 
for larger N we observe a range of values of k for which 
a stable packing cannot be found (a stable configuration 
corresponds to a small value of neighbor stability and 
hence a darker coloration in Fig. [2|c)). Nonetheless, as 
k increases, more stable configurations dominate. Con- 
versely, greater vorticity (indicated by a smaller numeri- 
cal value, and hence darker coloration, as this number is 
the cosine of the angular change) is evident for smaller N 
and larger k — and coincides with the values of neighbor 
stability which are stable. That is, a stable configuration 
will start to rotate, whereas an unstable one does not. 
Both travel and movement (as denned above) indicate 
similar results — stable configurations rotate and hence 
move more. There is no pattern in the ratio of travel to 
movement (not shown) as both quantities provide very 



similar information. 

In Fig. [3] we depict the phase transition more clearly 
by illustrating the value of these quantities as a function 
of only one of k or N. This figure displays the same 
data as in Fig. |2j when viewed together it is clear that 
one of two distinct behaviors persists in this system. Ei- 
ther a rigid fixed near neighbor structure with a rotating 
swarm and overall meander or, a swarm in which indi- 
viduals are constantly vying for position and lacking ro- 
tation or general drift. That is, we observe a transition 
from a rigid crystalline configuration (with each particle 
fixed in its relative position to its neighbors) to a gaseous 
phase (with particle constantly rearranging and moving 
past one another. The rigid configuration also exhibits a 
rotational force and hence a general meander. Nonethe- 
less, the gaseous phase is only transitional, by increasing 
k further (for fixed TV) the crystalline configuration re- 
sumes. 

While we have identified a transition to instability 
which can occur for a range of TV and k, the stable config- 
urations are both deterministic and rigid. In this respect 
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(a) local sep. (b) global sep. 
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FIG. 4. Followers For fixed k and N we now varying the 
number of followers k in the swarm. Here we illustrate with 
results for N = 50 and k = 6. We set r = 0, 1, 2, ... 50 and for 
each fixed r simulate, from random initial conditions a swarm 
for 20000 time steps, reporting dynamical measures over the 
lat 500 time step (sampling every 10-th). Mean and standard 
deviation of 30 trials are shown. 

they are, perhaps, not a good proxy for real animal be- 
havior. Hence, it is natural to ask what would need to be 
added to the model to introduce the natural variability 
of the real world. Stochastic perturbation is an obvious 
answer which we choose to ignore. Rather, we focus of 
whether introducing a small amount of aberrant (inho- 
mogeneous) behavior is sufficient to make a previously 
stable configuration unstable. While the effect of leaders 
in swarm simulations has been widely studied (particu- 
larly from the perspective of control) we choose to intro- 
duce aberrant "followers" . In a model such as the one we 
are looking at here, it is fairly straightforward to see that 



a sentient leader will easily be able to direct the motion 
of the flock. Rather, we ask what happens if we introduce 
a small number of rigid neighbor links. That is, rather 
that the neighbors of being the particles x^ such 
that dij = \\x^ — x^\\ is least, we simply prescribe that 
particle sc^M) is and always must be one of the neighbors 
of x^\ The distance da^ may not necessarily be small- 
est, particle is a "neighbor" of x^ only for the 
purposes of computing the update rules ([2J3]). Moreover, 
the relationship is not necessarily reciprocated: particle 
a^M) has no compulsion to consider x^ as a neighbor 
— unless they actually are close. Finally, we do this for 
a fixed small number of neighbors r <C Nk. If r = k then 
each particle is a follower of exactly one other random 
particle, iff r > k will a particle be a follower of more 
than one other particle. 

Figure [4] depicts the effect of introducing these follow- 
ers to our modeling scheme. Local separation shows a 
clear linear decrease with increasing r. Global separation, 
however, increases dramatically. Beyond about r = 12 
the swarm will fracture and form separate sub-groups. 
This behavior persists until r > 43. Most dramatically, 
from Fig. |4jc) we see a sudden transition at r — 10 
- from a rigid solid state to a gaseous one. Increasing 
r — >• 50 eventually reverses this trend. Finally, from Fig. 
|4](d-f) we see that as the gaseous phase occurs the swarm 
ceases to rotate and overall movement is less. 

We have shown that movement interaction rules 
adapted from nature can exhibit a wide range of inter- 
esting collective phenomena. The swarming behavior is 
notably coherent (forming a single flock) and exhibits a 
phase transition from a rigid crystalline state to a gaseous 
fluid-like one as a function of flock size and local neigh- 
borhood size. Nonetheless, the stable configuration is 
dominant and will re-emerge by further increasing the 
neighborhood size k. For a given stable configuration 
(fixed N and k) we demonstrate that the introduction 
of a small number of followers leads to a threshold phe- 
nomena. Beyond a moderate value of r (N = 50 and 
k = 6 we demonstrate that r = 10) the swarm behavior 
again transitions to a gaseous phase. This persists until 
for very large r the flock again becomes stable, possi- 
bly fracturing. Notably, the various dynamical behaviors 
we have observed here are deterministic. We have been 
able to offer a simple interaction model which mimics the 
variability found in nature without resorting to a proba- 
bilistic model. It would, perhaps, be interesting to extend 
this system to the case where the followers are not config- 
ured randomly but target specific individuals — perhaps 
following a scale-free distribution. 
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